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SECTION  1 


INTRODUCTION 

The  purpose  of  this  report  is  to  bring  to  light  an 
often  overlooked  physical  effect  that  may  be  a  factor  in 
the  evolution  of  electrostatic  plasmas.  This  effect  can 
alter  their  morphological  evolution  and  may  also  be  related 
to  if  not  responsible  for  their  structuring. 

Most  electrostatic  analyses  and  simulations  of 
plasma  clouds  assume  that  the  coupling  between  the  neutrals 
and  ions  is  basically  one-way.  They  account  for  the 
influence  of  neutral  motion  on  ion  motion  but  neglect  how 
the  reaction  forces  alter  the  flow  fields  of  the  neutrals. 

As  shown  in  Section  2  of  this  report,  indeed,  these 
reaction  forces  are  significant  and  can  alter  the  evolu¬ 
tion  of  the  plasma.  Further,  turbulence  associated  with 
these  forces  may  be  the  cause  of  ion  cloud  structuring. 
These  effects  are  of  direct  relevance  to  both  HANE  plasmas 
as  well  as  barium  plasma  clouds. 

In  Section  3  of  this  report  the  results  of  numerical 
simulations  which  attempt  to  model  the  influence  of  the 
ions  on  the  neutrals  as  well  as  the  influence  of  the 
neutrals  on  the  ions  is  presented.  Two  different  plasma 
cloud  cases  are  considered.  One  case  demonstrates  the 
effects  of  weak  coupling  parameters  similar  to  those 
believed  to  characterize  the  barium  cloud  release  Avefria 
Dos.  The  other  case  attempts  to  demonstrate  more 
significant  coupling  with  parameters  similar  to  the  ESTHER 
release  of  STRESS.  McDonald  et  al. ,( 1979)  notes  that  there 
are  distinct  large  scale  morphological  differences  between 
these  two  barium  clouds  and  presents  one  explanation  of 
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these  differences.  The  results  of  the  numerical 
simulations  discussed  in  Section  3  suggest  that  coupling 
between  ions  and  neutrals  may  be  an  additional  source  for 
these  differences. 

Conclusions  of  the  report  are  presented  in  Section 

4. 
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SECTION  2 


COUPLED  ION-NEUTRAL  FLOW  CONSIDERATIONS 


2 . 1  MOMENTUM  TRANSFER 


The  force  equation  for  ions  embedded  in  a  neutral 
gas  under  the  influence  of  electromagnetic  forces  is  given 
by 


m.N 

1 


.  3^  =  eN.  (E  +  u  X  B) 
1  dt  1 


where  u,  N.,  m.  and  u  are  the  ion  velocity,  particle 
density,  particle  mass,  and  particle  mobility  in  the 
neutral  gas;  where  e,  E,  and  B  are  the  charge  per 
particle,  electric  field  vector  and  magnetic  field 
vector,  respectively;  and  where  v  is  the  neutral  gas 
velocity.  The  second  term  on  the  right-hand  side  of  the 
equation  represents  the  forces  that  the  neutral  wind 
exerts  on  plasma  clouds  and  the  momentum  transfer  from 
the  neutral  gas  to  the  ions. 

Since  momentum  is  conserved  the  momentum  gained  by 
the  ions  must  represent  a  loss  in  momentum  by  the  neutrals. 
The  force  equation  for  the  neutrals  with  this  reaction 
term  is  thus 


dv  2 

mN  ^  =  -Vp  +  mNvV  v  + 


where  N,  m,  p,  and  v  are  the  neutral  density,  neutral 
particle  mass,  pressure,  and  kinematic  viscosity, 
respectively.  This  equation  is  the  Navier-Stokes 
equation  with  the  inclusion  of  the  ion  drag  force. 

Of  significant  interest  is  the  magnitude  of  the  ion 
drag  term.  The  absolute  value  of  the  drag  deceleration 
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can  be  obtained  by  dividing  the  neutral  force  equation  by 


mN  giving 


dv  Vp  .  „2 

dt  mN  mN  p 


(u  -  V) 


and  focussing  attention  on  the  last  term  on  the  right-hand 
side  of  the  equation.  This  term  is  essentially  the  product 
of  a  mass  density  ratio  between  the  ions  and  neutrals,  the 
slip  velocity,  and  the  collision  frequency  (defined  by 
e/ym^) . 

For  practical  purposes  it  is  easier  to  consider  the 
form  presented  above  because  e/mNy  can  be  treated  as 
roughly  constant  with  altitude.  From  McDaniel  (1964,  p. 
725)  the  ionic  mobility  is  given  by 


(k-1)  mN 


m  +  m A  1/2 


where  k  is  the  dielectric  constant  of  the  neutral  gas. 
The  term  k-1  should  be  closely  proportional  to  neutral 
particle  density  especially  at  ionospheric  particle 
density  levels  and  the  mobility  consequently  takes  the 
form 


~  N 


m  +  m.\l/2 


where  A  is  independent  of  density  and  particle  mass. 

The  dependence  of  y  on  1/N  makes  the  denominator  of  the 
ion  drag  term  independent  of  neutral  particle  density. 
While  y  and  N  are  strongly  dependent  upon  altitude 
the  only  altitude  dependence  of  the  drag  coefficient 
e/mNy  is  due  to  the  slowly  changing  effective  particle 
mass  and  particle  polarization  of  the  air  mixture. 


The  coefficient  e/mNy  for  barium  ions  can  be 

evaluated  from  published  laboratory  measurements  of  the 

-4  2 

mobility  of  2.25  x  lo  m  /volt-sec  for  bariiim  in  N2 
at  1  atm  at  20“C  (Brown,  1959,  p.  62) .  With  mN  as 
1.23  kgm/m^  (Knapp  and  Schwartz,  1975,  p.  8-2)  at  20 ®C 
at  1  atm 

e/itiNy  =  5.78  x  10  m^/sec 


The  significance  of  the  value  of  e/Nmy  can  be 
seen  by  considering  stationary  ions  and  by  neglecting  all 
other  neutral  wind  forces.  Integrating  the  Navier-Stokes 
equation  in  time  gives 


/ 


dv 

dt 


dt 


N. 

1 


vdt 


or 


ymN  ^Ti 


where  Av  is  the  change  in  velocity  of  a  neutral  fluid 
particle  passing  through  a  cloud  of  ions  with  integrated 
content  of  Thus  e/Nmy  is  the  change  in  neutral 

velocity  per  unit  of  integrated  plasma  ion  content.  For 
barium  the  neutral  wind  velocity  is  thus  decreased  by 
1  m/sec  for  every  1.73  x  10^^  ions  per  meter  squared  it 
passes  through  if  the  pressure  and  viscous  forces  are 
negligible. 

For  small  barium  clouds,  such  as  Avefria  Dos,  the 

peak  integrated  content  that  the  wind  passes  through  may 
16  2  7 

be  only  3  x  10  el/m  (10  el/cc  over  a  1  km  e-fold 
radius)  producing  only  a  10  m/sec  change  in  the  neutral 
wind  speed.  If  the  ion-neutral  slip  speed  is  high,  say 
60  m/sec,  a  10  m/sec  change  in  wind  speed  may  not  produce 
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a  noticeable  difference  in  evolution  from  the  case  of  no 

drag  effects.  For  larger  clouds  (48  kgm  at  180  km)  the 

peak  field  line  integrated  TEC  can  be  as  large  as 
17  2 

2.4  X  10  el/m  (Linson  and  Bcixter,  1977)  at  early  time. 
At  early  time  in  the  barium  cloud  development  these  large 
values  of  TEC  can  exist  in  the  dimension  transverse  to 
the  wind.  Changes  in  the  wind  velocity  as  large  as  130 
meters  per  second  are  predicted  (if  other  forces  are 
ignored)  which  should  significantly  affect  the  ion  cloud 
evolution. 

2.2  NEUTRAL  FLOW  CONSIDERATIONS 

The  problem  of  the  flow  of  ionospheric  neutrals 
through  a  cloud  of  ions  fixed  to  field  lines  is  analogous 
to  the  flow  of  a  gas  through  a  porous  medium  such  as  a 
sponge.  While  there  is  no  analytical  solution  to  this 
problem  the  analogy  proves  to  be  useful  when  considering 
the  properties  of  the  neutral  flow  to  be  expected. 

To  first  order  the  neutral  flow  is  reduced  in 
proportion  to  the  total  ion  content  passed  through.  This 
first  order  solution  would  result  from  a  balance  between 
the  fluid  momentum  term  and  the  drag  term  in  the  Navier- 
Stokes  equation  with  other  terms  ignored.  The  neutral 
flow  resulting  from  this  type  of  solution  however  is  not 
divergenceless.  Consequently,  the  density  is  enhanced 
and  reduced  in  different  regions  of  the  ion  cloud  and 
pressure  effects  start  to  influence  the  flow.  The 
pressure  effects  cause  the  flow  to  be  diverted  around  the 
ion  cloud  and  to  be  pulled  back  into  the  wake. 


A  significant  force  in  the  higher  altitude  neutral 
atmosphere  is  viscosity.  Viscosity  tends  to  smooth  the 
flow  pattern  resulting  from  the  pressure  and  drag  forces. 
While  the  flow  pattern  remains  similar  in  form,  the 
magnitude  of  the  decrease  in  the  neutral  velocity  in  the 
ion  cloud  wake  is  reduced  by  viscous  forces.  Viscous 
forces  thus  tend  to  de-emphasize  the  effects  that  ion 
drag  has  on  plasma  cloud  evolution,  especially  as  the 
cloud  altitude  is  increased. 

Perhaps  one  of  the  more  exciting  facets  of  ion 
neutral  coupling  effects  is  that  neutral  turbulence  may 
be  created.  One  would  expect  the  flow  through  a  porous 
obstacle  to  be  turbulent  especially  on  the  wake  side  of 
the  obstacle.  Note  that  this  side  of  the  ion  cloud  is 
the  side  that  typically  forms  structures  suggesting  that 
neutral  turbulence  could  be  the  source  of  ion  cloud 
structure.  Further  evidence  for  this  suggestion  is  that 
the  size  of  the  eddies  in  the  neutral  flow  roughly 
corresponds  to  the  sizes  of  the  ion  cloud  structures 
observed.  The  neutral  turbulence  structuring  mechanism 
offers  an  alternative  to  the  mechanism  discussed  by 
McDonald  et  al.  (1981)  to  explain  structure  growth  and 
size  stcibilization. 

The  inferences  about  the  possibility  of  turbulence 

and  the  size  of  their  eddies  are  made  from  Reynolds  number 

considerations.  Basically  the  size  at  which  the  viscosity 

terms  in  the  Navier-Stokes  equations  are  equal  to  the 

momentum  terms  of  the  equation  can  be  used  as  a  smallest 

eddy  size  (Landau  and  Lifshitz,  1959,  p.  122)  with  an 
1/4 

adjustment  of  R  '  to  the  slip  velocity  where  R  is 
the  Reynolds  number  for  the  entire  cloud.  The  momentum 
terms  can  be  approximated  to  be  of  magnitude 
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V  Av 
L 

where  L  is  the  size  of  the  structure,  v  is  the  slip 
speed  through  it,  and  v  is  the  change  in  slip  speed. 
The  viscosity  term  can  be  similarly  estimated  to  be 


where  v  is  the  kinematic  viscosity.  The  two  terms 
balance  when  the  Reynolds  number  for  the  flow  through 
the  feature  is  1,  i.e.,  when. 


By  assuming  an  arbitrary  ion  neutral  slip  speed  of  100 
meters  per  second  through  the  entire  cloud  the  eddy  size 
can  be  given  in  terms  of  the  kinematic  viscosity.  The 
kinematic  viscosity  is  a  strong  function  of  altitude  and, 
consequently,  so  is  the  eddy  size.  The  eddy  size  for  a 

range  of  altitudes  is  presented  in  Table  1. 

1/4 

The  R  '  correction  to  the  eddy  size  has  been 
made  with  the  assumption  that  the  overall  plasma  cloud 
size  is  5  kilometers.  Essentially  this  correction  term 
arises  because  the  effective  flow  speed  differential  Av, 
through  the  eddy  region  is  somewhat  less  than  the  expected 
flow  speed  differential  through  the  total  cloud. 

The  proposed  mechanism  for  creating  structure  has 
many  desirable  features  in  that  it  does  not  rely  upon 
seed  structure  in  the  ionosphere  to  produce  growth  as 
does  the  gradient  drift  (Linson  and  Workman,  1970) 
instability.  It  does  have  its  drawbacks,  however,  which 
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should  be  pointed  out.  One  drawback  is  that  the  Reynolds 
number  required  to  establish  turbulence  may  be  too  large. 
Landau  and  Lifshitz  (1959,  p.  105)  mention  that  for  the 
case  of  a  solid  cylinder  the  Reynolds  number  required  for 
the  flow  to  become  turbulent  is  34.  For  the  case  of  the 
ion  drag  effect  perhaps  such  a  large  Reynolds  number  is 
not  required  because  of  the  porous  nature  of  the  flow  but 
in  any  event  the  required  Reynolds  number  should  probably 
be  greater  than  ten.  Apparently,  for  clouds  less  than  5 
kilometers  in  extent  the  ions  must  be  at  an  altitude  less 
than  160  kilometers  to  produce  turbulence. 

Another  drawback  to  note  is  that  the  eddy  sizes 
predicted,  though  close  to  the  sizes  of  barium  structures, 
are  nevertheless  somewhat  large  especially  for  the  higher 
altitude  clouds.  For  many  clouds  the  slip  velocity  is 
closer  to  25  meters  per  second  in  which  case  the  uncorrected 
eddy  sizes  increase  by  a  factor  of  4. 

The  high  kinematic  viscosity  in  the  180  to  200 
kilometer  altitude  regime  makes  structuring  by  neutral 
turbulence  to  the  scale  sizes  observed  in  barium  releases 
difficult.  It  may  be  the  case,  however,  that  electrical 
coupling  to  lower  (130  kilometer)  altitude  turbulence 
could  play  a  role  in  setting  structure  sizes.  This  thought 
leads  to  still  another  structuring  mechanism.  Perhaps 
barium  structure  is  created  by  neutral  turbulence,  either 
pre-existing  or  stimulated  by  images,  in  the  130  kilometer 
altitude  regime  which  is  coupled  electrostatically  to  the 
barium  cloud  along  field  lines.  The  100  meter  size  eddies 
predicted  in  this  regime  from  viscosity  considerations  are 
in  good  agreement  with  observed  striation  sizes. 

In  summary,  the  neutral  flow  through  a  plasma  cloud 
is  similar  to  the  flow  of  a  gas  through  a  porous  medium. 
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Because  the  Reynolds  number  for  the  flow  can  be  large, 
the  neutral  flow  produced  may  be  turbulent.  Neutral 
turbulence  is  more  likely  at  lower  altitudes  where  the 
kinematic  viscosity  is  less  and  probably  difficult  to 
produce  above  200  kilometers.  Turbulence  in  the  neutral 
flow  could  conceivably  be  an  alternative  source  of 
electrostatic  plasma  cloud  structuring. 


2.3  FIRST  ORDER  INFLUENCE  OF  ION  DRAG  EFFECTS  ON  PLASMA 

CLOUD  MOTION 

Neutral  flow  through  ionospheric  plasma  creates  an 
electrical  current,  J  ,  carried  by  ions  with  value 
specified  by 

eN.  (uB)  eN. 

J  =  — - - -  V  X  i  +  - i - y  V 

1  +  (yB)^  ^  1  +  (yB)^ 

If  the  current  is  divergenceless  it  will  not  produce 
motion  in  the  plasma  other  than  that  required  of  the  ions 
to  carry  the  current.  On  the  other  hand  if  the  current 
is  not  divergenceless  electrical  fields  will  be  created 
to  bleed  off  charge  build-up.  These  electric  fields 
become  the  source  of  plasma  motion. 

Divergent  wind  currents  in  plasma  clouds  are 

typically  created  by  variations  in  the  plasma  density. 

Divergent  wind  currents  can  also  be  produced  by  variations 

in  the  neutral  wind  transverse  to  the  geomagnetic  field 

even  though  the  plasma  density  is  constant.  While  the 

divergence  of  the  neutral  wind  is  typically  close  to  zero, 

the  neutral  wind  current  depends  most  strongly  upon 

V  •  V  X i  which  is  typically  not  zero  in  ion  drag  flows, 
z 
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Ion  clouds  which  produce  deceleration  of  the 
neutral  flow  can  stimulate  both  of  these  two  different 
sources  of  divergent  wind  currents — a  situation  which 
leads  to  interesting  plasma  flow  properties.  Consider 
the  neutral  flow  through  an  ion  cloud  with  a  Gaussian 
spatial  density  distribution.  To  first  order  the  neutral 
flow  in  the  wake  is  linearly  related  to  the  integrated 
ion  content  that  the  neutrals  pass  through.  Neglecting 
pressure  and  viscosity  effects  the  idealized  flow  is  as 
shown  in  Figure  2-1.  The  flow  in  this  figure  is  always 
directed  toward  the  right,  only  its  magnitude  in  this 
direction  changes.  In  the  interior  region  of  the  plasma 
cloud  the  flow  is  not  divergenceless  and  hence  its  pattern 
is  somewhat  unrealistic.  In  the  wake  region  however  the 
flow  is  divergenceless  and,  thus,  potentially  realistic. 

The  divergence  of  the  wind  driven  component  of  the 
current  is  proportional  to 


for  the  case  where  the  magnetic  field  is  into  the  paper 
and  where  jjB  has  been  assumed  to  be  large.  The  regions 
of  charge  accumulation  for  this  situation  are  schematically 
shown  in  Figure  2-2.  In  addition  to  the  customary  charge 
accumulation  at  the  top  and  bottom  of  the  ion  cloud  is  a 
charge  accumulation  associated  with  the  wake  of  the 
cloud. 

The  plasma  flow  associated  with  the  customary  charge 
accumulation  is  such  that  the  ion  cloud  is  moved  toward 
the  right  along  the  centerline  of  the  cloud.  The  flow 
associated  with  the  wake,  on  the  other  hand,  is  directed 
toward  the  left  along  the  wake  centerline.  Background 


Figure  2-1.  Idealized  neutral  flow  through  ion 
cloud.  Dotted  lines  represent  iso¬ 


density  contours  of  the  ion  cloud  and 
solid  lines  represent  iso-speed  contours 
of  the  x-directed  neutral  wind. 


Figure  2-2.  Regions  of  charge  accumulation  created 

by  divergences  in  the  wind  current  shown 
superimposed  on  flow  and  density  contours 


17 


plasma  in  the  wake  is,  thus,  directed  toward  the  center 
of  the  barium  cloud.  This  flow  is  somewhat  unusual  in 
that  plasma  is  moving  toward  the  left  in  response  to 
winds  toward  the  right.  This  phenomenon  occurs  bascially 
because  the  plasma  is  trying  to  maintain  a  constant  ion- 
neutral  slip  speed  in  order  to  maintain  a  constant 
y-directed  current. 

One  consequence  of  the  mixture  of  the  two  different 
types  of  flow  is  that  a  dimple  at  the  front  of  the  ion 
cloud  can  be  produced  if  the  counter  directed  wake  flow 
is  strong  enough.  The  plasma  motion  toward  the  right  at 
the  top  and  bottom  of  the  cloud  may  become  faster  than  at 
the  centerline  and  as  a  result  the  cloud  may  be  deformed 
until  it  begins  to  bifurcate.  This  type  of  behavior  has 
not  been  reproduced  in  simulations  which  do  not  include 
ion-neutral  coupling  and  is  certainly  of  interest  for 
understanding  structuring  processes  in  plasma  clouds. 

Its  effect  is  apparent  in  one  of  the  simulations  which 
follow. 


SECTION  3 


COUPLED  ION-NEUTRAL  FLOW  NUMERICAL  SIMULATION 


3.1  COUPLED  ION-NEUTRAL  FLOW  RELATIONS 


In  the  following  section  results  of  simulations  of 
the  evolution  of  plasma  clouds  which  include  ion-neutral 
coupling  effects  are  presented.  In  the  simulations  the 
ion  drag  effects  on  the  neutral  flow  are  included.  This 
section  presents  the  relations,  the  assumptions,  and  the 
approximations  used  for  the  simulation  calculations. 

Five  relations  are  used  in  the  calculations  which 
involve  five  different  quantities,  namely,  the  neutral 
density,  N,  the  plasma  density,  N.,  the  flow  potential, 
0,  and  the  x-  and  y-components  of  the  neutral  wind  v. 

The  relations  solved  to  find  the  ion  behavior  are 
the  same  as  those  conventionally  used  in  2-dimensional 
electrostatic  plasma  evolution  simulations.  A  split 
step  scheme  is  used.  The  plasma  flow  field  is  solved  in 
the  potential  solving  step  from  the  specified  values  of 
plasma  density  and  neutral  wind.  The  flow  is  then  used 
to  update  the  plasma  density  values  in  the  plasma  convec¬ 
tion  step.  The  plasma  convection  step  numerically 
advances  the  plasma  density  by  using  the  continuity 
relation 


N.  V(b  X  i 
i  z 


Note  that  the  plasma  flow  velocity  u  is  assumed  to  be 
incompressible  and  <P  is  the  potential  function 
describing  this  flow; 

u  =  V(J)  X  i 

z 
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The  potential  solving  step  solves  the  potential  equation 
which  specifies  that  the  divergence  of  the  total  current 
be  zero.  For  the  simulations  performed  the  Hall  terms  are 
assximed  to  be  negligible.  This  assumption  is  frequently 
called  the  F-layer  assumption  and  it  implies  the  following 
relation  for  the  potential: 


V*N.  V*  =  -V'N.v^i 
1  1  z 

In  this  relation  the  neutral  wind  N.  and  v  are  assumed 

1 

known  and  <t)  is  solved  for  using  matrix  inversion  tech¬ 
niques  . 

Once  and  <t)  are  known  the  effect  of  ion  drag 

on  the  neutral  flow  can  be  calculated.  The  neutral  density 
and  flow  are  also  advanced  through  time  using  a  split  step 
technique.  Given  a  known  neutral  flow  velocity  the  neutral 
density  is  advanced  in  time  using  the  continuity  equation 


9N 

It  = 


-V 


Nv 


The  neutral  density  flow  is  updated  through  time  using  the 
Navier-Stokes  equation  with  an  isothermal  ideal  gas 
approximation  to  the  pressure  term  giving 


=  -(vV)v  -  —  Vin  N  +  vV^v  +  a(N.  -  N,  )  (u  -  v) 
dt  m  lb 

where  Nj^  is  the  background  ion  density,  k  is  Boltzmann's 
constant,  T  is  the  neutral  gas  ambient  temperature,  m  is 
the  neutral  particle  mass,  v  is  the  kinematic  viscosity 
and  a  is  the  coupling  coefficient  e/mNy  discussed  in 
Section  2.1.  The  above  equation  constitutes  two  relations, 
one  relation  for  each  component  of  v. 
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while  it  is  probably  true  that  the  adiabatic 
assumption  for  the  slightly  compressible  component  of 
the  neutral  flow  is  probably  closer  to  the  actual 
physical  situation,  the  isothermal  assumption  provides 
a  very  accurate  approximation  because  the  deflections 
of  the  density  from  its  aunbient  value  are  very  small 
(less  than  1%) .  Use  of  the  isothermal  approximation 
leads  to  terms  which  are  more  straightforward  to  evaluate 
in  the  Navier-Stokes  equation.  If  the  adiabatic  assump¬ 
tion  were  made  then  terms  involving  the  evaluation  of 
the  density  raised  to  a  non-integral  power  would  require 
evaluation. 

Subtraction  of  the  background  density  from  the  drag 
component  is  necessary  to  remove  the  influence  of  drag 
effects  on  the  background.  The  tendency  for  the  ambient 
neutral  wind  to  decelerate  in  the  ionosphere  is  probably 
balanced  by  pressure  gradients  and  thus  a  removal  of  the 
background  ion-neutral  drag  seems  appropriate. 

In  order  to  assure  accuracy  and  stability  (Ames, 
1977,  p.  264)  the  time  step  size  used  for  the  neutral 
flow  calculation  must  be  fairly  small.  The  time  step 
must  be  so  small  that  flow  at  the  fastest  expected 
neutral  wind  velocity  plus  the  sound  speed  will  not  cross 
more  than  70%  of  a  grid  cell.  In  contrast  the  time  step 
used  for  the  ion  flow  calculation  can  be  significantly 
larger.  The  time  step  must  only  be  small  enough  so  that 
flow  at  the  fastest  ion  speed  will  not  cross  more  than 
50%  of  a  grid  cell.  Typically,  in  a  calculation  in  the 
frcune  of  reference  of  the  background  ions  the  fastest 
ion  speed  will  be  less  than  the  neutral  speed.  The 
disparity  between  the  time  step  required  by  the  neutral 
and  ion  flow  motivates  the  use  of  two  different  time 
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step  values.  In  the  calculational  results  to  be  presented 
the  neutral  calculation  is  iterated  numerous  (16  to  50) 
times  for  every  iteration  of  the  ion  flow. 

Finite  difference  techniques  were  applied  using  a 
rectangular  grid  of  square  grid  cells.  The  problems  were 
formulated  to  be  symmetric  around  the  horizontal  axis  and 
advantage  was  taken  of  this  symmetry  to  reduce  storage 
requirements.  The  sampled  functions  are  assumed  tc  have 
a  planar  dependence  over  each  cell. 

The  boundary  values  of  ion  and  neutral  density  and 
the  X-  and  y-components  of  the  neutral  wind  are  all 
assumed  to  be  such  that  the  quantities  have  continuous 
normal  derivatives  across  the  boundary.  The  x-component 
of  the  neutral  wind  is  further  assumed  to  have  a  fixed 
value  on  the  left  hand  boundary  that  is  representative  of 
the  incoming  neutral  flow  speed.  The  boundary  value  of 
the  ion  flow  potential  is  set  such  that  the  ion  motion 
interior  to  the  boundary  is  not  influenced  by  the 
boundary's  presence.  The  potential  equation  is  effectively 
solved  with  the  free  space  Green's  function  for  the  problem. 

The  initial  conditions  for  the  neutral  parameters 
are  that  the  neutral  density  and  wind  components  are 
constant.  The  x-component  of  the  wind  is  assumed  to  take 
the  value  of  the  ion-neutral  slip  speed.  The  y-component 
of  the  wind  is  assumed  to  be  initially  zero.  The  neutral 
density  is  assumed  to  have  an  initial  value  of  unity. 

The  initial  value  of  the  ion  density  is  assumed  to 
have  a  Gaussian  spatial  variation  superimposed  on  a 
constant  background  of  unity.  The  peak  value  of  the 
Gaussian  is  three  in  both  cases,  two  larger  than  the 
assumed  background  ion  density.  This  value  causes  the 
ion-neutral  slip  to  be  roughly  one-half  of  the  wind 
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velocity  which  is  in  fair  agreement  with  observations  of 
these  clouds.  The  initial  e-fold  radius  has  been  set 
according  to  which  release  is  modelled  as  specified  in 
Table  3-1.  Table  3-1  also  lists  many  other  important 
parameters  used  in  the  two  simulations,  the  results  of 
which  are  discussed  in  the  following  section. 


3.2  NUMERICAL  RESULTS  FOR  THE  SIMULATIONS  OF  AVEFRIA 

DOS  AND  IRIS 

The  results  of  a  coupled  ion-neutral  interactive 
flow  calculation  for  a  small  (1.7  km  e-fold)  dense 
(3.5  X  10^  el/cc  peak)  plasma  cloud  are  illustrated  in 
Figures  3-1  thru  3-13.  The  parameters  used  in  the 
calculation  are  presented  in  Table  3-1  and  were  chosen 
to  characterize  the  pressure  and  viscosity  at  a  200 
kilometer  altitude.  Above  each  figure  the  contour 
values  are  given  and  the  contours  are  coded  such  that 
the  highest  value  contour  is  bold  and  the  lowest  value 
contour  is  dotted.  The  distance  values  in  the  labels 
are  assumed  to  be  in  kilometer  spatial  units. 

Figure  3-1  shows  contours  of  the  initial  plasma 
density.  The  initial  density  profile  has  a  Gaussian 
spatial  dependence  which  is  centered  to  the  left  hand 
side  of  the  grid.  The  neutral  wind  is  initially  directed 
toward  the  right  and  centering  the  cloud  to  the  left 
allows  for  advancement  of  the  simulation  to  later  times 
without  moving  the  frame  of  reference.  The  neutral  flow 
is  calculated  for  the  first  sixteen  neutral  time  steps 
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Table  3-1 


Simulation  Parameters 


Avefria  Dos 

IRIS 

Units 

Initial  Values 

Radius — R 

o 

1.7 

5 

km 

Neutral  wind — v 

X 

125 

55 

m/sec 

Neutral  Parameters 
and  Coefficients 

Pressure — kT/m 

3.4  X 10^ 

2.4  xlO^ 

2  2 
m  /sec'^ 

Viscosity — v 

1.6  X 10^ 

5.3  X  lo'^ 

2  , 

m  /sec 

Coupling — a 

1  X lo"^^ 

-12 

6  xio 

3  , 

m  /sec 

Peak  density 
(implied  by  a) 

3.5  X 10^ 

2  xlO^ 

particles/cc 

Reynolds  number 
(R^v/v) 

1.3 

5.2 

Sampling  Parameters 

Grid 

42  X  24 

32  X  28 

X  by  y 

Total  Distance 

20  X 11.4 

47  X  41.1 

km  in  X  by 
km  in  y 

Ax  =  Ay 

476 

1,469 

meters 

At . 
ion 

2 

10 

seconds 

At.  /At  ^  , 

ion  neutral 
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Figure  3-1.  Initial  plasma  density  for  Avefria  Dos 
calculation. 
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Figure  3-2.  The  x-component  of  the  neutral  wind  after 

10  seconds  of  initial  evolution.  Units  are 
km/second. 
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Figure  3-3.  The  y-component  of  the  neutral  wind  after 
10  seconds  of  initial  evolution. 
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Figure  3-4. 
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The  neutral  density  after  10  seconds  of 
initial  evolution. 
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Figure  3-5.  The  electrostatic  potential  resulting  from 
initial  plasma  density  and  initial  wind 
values. 


Figure  3-6.  The  x-component  of  the  initial  ion  velocity. 
Units  are  kiti/sec. 
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Figure  3-7.  The  y-component  of  the  initial  ion  velocity. 
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Figure  3-8. 
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The  plasma  density  after  150  seconds  of 
evolution. 
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Figure  3-9.  The  x-component  of  the  neutral  wind  after  150 
seconds  of  evolution.  The  wind  drag  effect 
creates  a  10  meter/second  decrease. 
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Figure  3-10.  The  y-component  of  the  neutral  wind  after 
150  seconds  of  evolution. 
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Figure  3-11.  The  neutral  density  after  150  seconds  of 
evolution. 


30 


l«  OCT  la  TINE  lEl 


Ml  HIM  -«.tf2n*oa  H*i  ».ifa«E-ea 

COMTOUN  -.<3IC-OI  -.lOK-OI  -.IME-OI  -.il5E-Sa 

•iiic-oa  .im-oi  .3o*c-oi  .oie-oi  .sj^e-oi 


Figure  3-12.  The  electrostatic  potential  after  150 
seconds  of  evolution. 
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Figure  3-13.  The  ion  density  after  150  seconds  of 
evolution. 
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assuming  the  plasma  velocity  is  zero.  Plots  of  the 
resultant  neutral  flow  parcuneters  are  shown  in  Figures 
3-2  thru  3-4.  The  neutral  wind  values  are  then  used  to 
determine  the  electrostatic  potential  shown  in  Figure 
3-5.  Because  the  free  space  boundary  condition  algorithm 
has  been  used  the  resultant  potential  is  not  affected  by 
the  fairly  close  grid  boundary.  The  plasma  velocity 
components  implied  by  the  potential  are  shown  in 
Figures  3-6  and  3-7.  The  maximum  value  of  the  x-component 
of  the  plasma  velocity  is  53  meters  per  second,  roughly 
one-half  of  the  initial  125  meter  per  second  neutral  wind 
speed.  The  ion  neutral  slip  speed  is,  thus,  roughly  70 
meters  per  second  through  the  maximum  ion  density  with 
higher  slip  speeds  occurring  elsewhere.  For  subsequent 
neutral  flow  iterations  the  ion  flow  velocities  are  used 
to  determine  the  drag  force,  instead  of  assuming 
stationary  ions  as  is  done  in  the  initial  time  step. 

Figures  3-8  thru  3-12  depict  the  calculated  ion 
and  neutral  parcuneters  after  150  seconds  of  evolution. 

The  evolution  depicted  shows  little  influence  of  the 
ion  neutral  coupling  effect  because  the  change  in  the 
neutral  wind  is  fairly  small.  The  calculation  indicates 
a  change  of  only  11  meters  per  second.  This  small  value 
disagrees  with  a  60  meter  per  second  value  predicted 
using  the  drag  effect  relations  of  Section  2.1,  however, 
those  relations  neglect  viscous  forces.  A  comparison  of 
viscous  forces  to  drag  forces  reveals  that  for  this 
problem  viscosity  is  the  velocity  determining  force,  as 
might  be  expected  because  of  the  low  Reynolds  nximber. 
Because  of  the  viscous  forces  Avefria  Dos  should  evolve 
in  a  manner  that  is  insensitive  to  ion-neutral  drag 


33 


effects  and  the  results  of  past  electrostatic  analyses 
including  McDonald  et  al.  (1979)  should  apply. 

It  may  be  worthwhile  to  point  out  a  discrepancy 
between  the  calculated  results  and  the  observations  and 
data  presented  in  McDonald  et  al.  (1979) .  In  the 
calculated  results  shown  in  Figure  3-8  an  elongation  of 
the  cloud  can  be  noted  that  is  of  order  2  1/2  to  1. 

Cloud  elongation  is  a  function  of  the  slip  velocity 
which  is  itself  a  function  of  the  peak-to-background 
ratio  and  the  wind  speed.  The  parameters  were  chosen 
according  to  McDonald  et  al.'s  (1979)  range  of  acceptable 
values  and  yet  the  cloud  elongation  is  much  more  severe 
than  that  observable  (=1  1/2  to  1)  in  the  optical  data 
also  published  in  McDonald  et  al.  (1979) .  The  implica¬ 
tion  of  the  discrepancy  is  that  the  peak-to-background 
ratios  found  acceptable  by  McDonald  et  al.  (1979)  based 
on  cloud  structuring  properties  may  not  be  reflective  of 
the  observed  cloud  dynamics  which  were  not  monitored  in 
McDonald  et  al.  (1979) .  Indeed,  peak-to-background 
ratios  of  2  or  less  (in  place  of  the  value  of  3  currently 
modelled)  would  be  appropriate  to  model  cloud  dynamics. 
These  low  values  lie  outside  the  acceptable  range 
presented  in  McDonald  et  al.  (1979) ,  suggesting  that 
the  structuring  of  Avefria  Dos  may  not  be  fully  explained. 

Figure  3-13  shows  the  plasma  density  after  250 
seconds  of  evolution.  By  this  time  nvimerical  difficulties 
are  becoming  apparent  however  there  is  still  no  evidence 
of  the  dimpling  effects  expected  from  ion-neutral 
coupling. 

The  results  of  a  coupled  ion-neutral  interactive 
flow  calculation  for  a  large  (5  kilometer  e-fold)  plasma 
cloud  at  a  lower  altitude  (170  kilometers)  are  shown  in 
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Figures  3-14  thru  3-34.  The  calculation  assumes  that  the 
peak  ion  density  is  2  x  10^  el/cc.  Note  that  at  early  time 
this  value  is  probably  appropriate  for  a  barium  cloud 
similar  to  the  IRIS  event  of  PLACES  (Linson  and  Baxter, 
1977)  (which  was  similar  in  many  respects  to  ESTHER  of 
STRESS  and  SPRUCE  of  SECEDE-II) ,  however,  at  later  time 
this  peak  value  is  decreased  through  diffusion  along 
the  geomagnetic  field  line — an  effect  not  modelled  in 
this  calculation. 

The  initial  values  of  ion  density  are  shown  in 
Figures  3-14.  The  values  of  the  neutral  wind  x-  and 
y-components  and  the  neutral  density  perturbation 
produced  by  10  seconds  of  flow  through  this  stationary 
plasma  distribution  are  shown  in  Figures  3-15  thru  3-17. 
Note  that  the  x-component  of  the  neutral  wind  drops  by 
5  meters  per  second  during  this  initial  interval.  The 
electrostatic  potential  and  the  x-  and  y-components  of 
the  plasma  velocity  are  shown  in  Figures  3-18  thru  3-20. 

The  qualitative  nature  of  these  contour  plots  is  very 
similar  to  the  initial  values  obtained  for  the  Avefria 
Dos  calculation  indicating  that  the  effects  of  ion- 
neutral  wind  drag  make  little  difference  to  the  initial 
cloud  evolution.  Note  that  the  initial  x-component  of 
the  plasma  velocity  is  21  meters  per  second  indicating 
a  slip  speed  of  roughly  35  meters  per  second  at  the 
peak  density  value.  This  slip  speed  is  in  agreement 
with  the  slip  speed  observed  during  IRIS. 

Figure  3-21  shows  the  plasma  density  after  150 
seconds  of  evolution.  At  this  time  only  a  small  amount 
of  contour  distortion  has  occurred.  Because  of  the 
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Figure  3-14.  Initial  ion  density  for  coupled  ion-neutral 
flow  calculation  of  IRIS. 
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Figure  3-16.  The  y-component  of  the  neutral  wind 
after  first  time  step. 
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Figure  3-18.  The  results  of  the  first  time  step 

electrostatic  potential  calculation. 
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Figure  3-24.  The  neutral  density  after  150  seconds 
of  evolution. 
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Figure  3-30,  Electrostatic  potential  at  155  seconds. 
Quadripole  nature  is  evident. 
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Figure  3-31.  Electrostatic  potential  at  245  seconds. 
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Figure  3-33.  Electrostatic  potential  at  425  seconds. 

The  dipole  associated  with  the  wake  has 
completely  dominated  the  dipole 
polarization  associated  with  the  cloud. 
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lesser  wind  velocity  this  cloud  evolves  slower  than 
Avefria  Dos. 


The  neutral  flow  parameters  at  150  seconds  are 
shown  in  Figures  3-23  thru  3-24.  These  quantities  are 
also  still  evolving.  The  decrease  in  the  x-component 
of  the  neutral  wind  is  only  9  meters  per  second  at  this 
time. 

The  plasma  density  and  x-component  of  the  neutral 
wind  at  330  seconds  are  shown  in  Figures  3-25  and  3-26. 
This  time  is  roughly  equivalent  to  the  150  second  time 
in  the  Avefria  Dos  evolution  based  upon  proportional 
scaling  to  the  neutral  wind.  The  distortion  of  the 
ion  cloud  is  now  significantly  different  from  that  of 
Avefria  Dos  because  of  the  ion-neutral  coupling  effect. 

At  the  330  second  time  the  neutral  wind  has  nearly 
reached  a  steady  flow  situation.  The  decrease  in  the 
neutral  wind  is  23  meters  per  second.  This  value 
reflects  a  balance  between  ion  drag  and  neuti  il  wind 
momentum  loss  and  to  a  lesser  extent  the  viscosity. 

In  Section  2.1  the  decrease  in  the  neutral  wind 
considering  only  the  momentum  term  is  given.  Simple 
approximations  to  the  flow  pattern  can  be  made  to 
derive  another  useful  relation  which  also  accounts  for 
viscosity  and  ion-neutral  slip.  The  resulting  relation 
is 

aN  j 

Av  =  - 5 -  V  ,,  . 

aN,  +  2v/R  +  V/2R 
d  o  o 

where  a  is  the  coupling  coefficient,  is  the  peak 

ion  density  minus  the  background  density,  R^  is  the 
e-fold  radius,  v  is  the  ambient  wind  speed  in  the 


freime  of  the  background  ions,  and  the  ion- 

neutral  slip  value  at  the  peak  density  position.  The 
factor  of  2  which  is  present  in  the  viscosity  term  of 
the  denominator  arises  from  the  second  derivative  of 
an  assumed  Gaussian  spatial  variation  in  the  y-direction 
for  the  x-component  of  the  neutral  wind.  The  factor  of 
2  which  appears  in  the  momentum  term  of  the  denominator 
accounts  for  the  fact  that  the  change  in  speed  occurs 
roughly  across  the  complete  diameter  of  the  cloud.  For 
the  results  of  Figure  3-26  the  Av  predicted  is  19  meters 
per  second  using  the  paraimeters  ct  =  .006,  =  2 , 

V  =  .053,  R  =5,  V  =  .055  and  v  , .  =35  meters  per 

second.  This  prediction  is  in  reasonable  agreement  with 
the  Av  of  23  meters  per  second  observed. 

Continued  evolution  of  the  ion  cloud  is  shown  in 
Figures  3-27  thru  3-29.  Evidence  of  the  dimpling  process 
becomes  apparent  at  390  seconds  as  shown  in  Figure  3-27. 
This  dimpling  becomes  more  pronounced  at  420  seconds  in 
Figure  3-28  and  remains  until  later  time.  Figure  3-29 
shows  the  predicted  ion  cloud  at  the  final  time  step 
calculated,  naimely,  900  seconds.  Note  that  numerical 
difficulties  have  arisen  as  evidenced  by  the  presence  of 
plasma  depletions  and  by  the  numerous  isolated  regions 
of  peak  density;  however,  no  significant  change  in  the 
nature  of  the  evolution  has  arisen.  Also  note  that  it 
is  not  realistic  to  expect  the  high  levels  of  ion- 
neutral  coupling  being  modelled  to  actually  be  present 
in  IRIS,  ESTHER,  or  SPRUCE  to  these  late  times  because 
diffusion  along  the  geomagnetic  field  decreases  the  peak 
density. 

Another  interesting  aspect  of  this  calculation  is 
the  evolution  of  the  electrostatic  potential.  Figure 
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3-18  shows  the  initial  potential  as  having  a  dipole-like 
nature.  The  wind  drag  effect  changes  the  potential 
behavior  into  a  quadripole  field  because  of  the  accumula¬ 
tion  of  wind  drag  related  charge  in  the  wake  of  the 
cloud.  Figures  3-30  thru  3-34  illustrate  the  evolution 
of  this  quadripole  field  with  time.  Note  that  the 
magnitude  of  the  peak  potential  associated  with  the  wind 
drag  wake  grows  continuously  until  by  425  seconds  the 
wake  potential  has  completely  dominated  the  polarization 
potential  of  the  ion  cloud  itself.  This  behavior  agrees 
in  principle  with  the  mechanisms  outlined  in  Section  2 . 3 
although  the  magnitude  of  the  wake  potential  may  be 
somewhat  surprising. 
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SECTION  4 


CONCLUSIONS 


In  Section  2  it  is  shown  through  momentum  conserva¬ 
tion  arguments  that  the  reaction  forces  of  ion  clouds  on 
the  neutrals  may  be  significant.  This  fact  motivates  the 
coupled  ion-neutral  flow  calculations  presented  in  Section 
3.  The  results  of  the  numerical  calculation  suggest  that 
ion-neutral  coupling  effects  may  be  a  significant  factor 
in  the  evolution  of  electrostatic  plasma  clouds.  They 
suggest  that  morphological  differences  between  small  barium 
clouds,  such  as  Avefria  Dos,  and  large  barium  clouds,  such 
as  IRIS,  ESTHER,  and  SPRUCE,  could  result  from  ion-neutral 
coupling  effects. 

In  Section  2  a  simple  formula  for  the  change  in 
neutral  wind  velocity  based  upon  a  balance  between  momentum 
terms  and  ion  drag  terras  is  derived  which  indicates  that 
ion  drag  effects  could  be  significant  in  both  large  and 
small  barium  ion  clouds.  Section  2  also  notes  that 
viscosity  can  be  an  important  parameter  in  the  flow  of 
the  neutrals  especially  as  the  altitude  increases  because 
of  the  increases  in  kinematic  viscosity.  This  observation 
is  borne  out  by  the  numerical  results  of  Section  3  which 
show  that  the  neutral  wind  drag  effect  for  the  Avefria 
Dos  cloud  is  controlled  by  viscous  forces.  Viscosity 
limits  the  decrease  in  neutral  wind  speed  to  11  meters 
per  second  for  this  cloud  where  a  60  meter  per  second 
speed  decrease  might  be  expected  from  momentum  considera¬ 
tions  alone.  An  expression  which  accounts  for  both 
viscous  and  momentum  terms  to  approximate  the  decrease 
in  the  neutral  wind  is  presented  in  Section  3 . 2  and  the 
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results  of  both  the  IRIS  and  the  Avefria  Dos  flow 
calculations  are  found  to  be  in  reasonable  agreement. 

The  latter  relation  and  the  calculated  results  emphasize 
the  importance  of  viscosity  in  high  altitude  neutral  flow 
considerations . 

The  morphological  difference  which  arises  because 
of  ion-neutral  coupling  effects  is  an  initial  increase  in 
the  ion  cloud  dimension  perpendicular  to  the  wind  on  the 
steepening  side  and  then  an  eventual  dimpling  of  this 
side.  The  dimpling  is  significant  because  it  is  the 
first  such  change  in  curvature  of  a  contour  noted  in 
well-resolved  unseeded  numerical  simulations. 

Perhaps  this  type  of  coupling  phenomenon  might  be 
responsible  for  the  two-pronged  fork  appearance  of 
PLACES  event  HOPE;  although,  it  should  be  noted  that  the 
morphology  resulting  from  ion-neutral  coupling  effects, 
as  shown  in  Section  3.2  Figure  3-29  in  particular,  differs 
from  the  two-prong  fork  morphology  of  HOPE  in  that  the 
angle  of  bifurcation  is  obtuse  in  the  numerical  results 
but  acute  in  the  HOPE  data.  Perhaps  better  agreement 
could  be  obtained  from  a  three-dimensional  analysis  of 
the  problem. 

In  Section  2 . 2  turbulence  in  the  neutral  wind  is 
suggested  as  a  possible  source  for  structure.  The  scale 
size  of  expected  eddy  cells  at  lower  altitudes  are  in 
close  agreement  with  the  structure  sizes  actually 
observed.  Since  few  other  physical  mechanisms  suggest 
the  proper  scale  sizes  of  structure  the  neutral  turbu¬ 
lence  mechanism  is  viewed  as  a  promising  explanation. 
However,  in  light  of  the  high  viscosity  of  the  higher 
altitude  regions  and  its  control  on  the  neutral  flow  as 
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as  observed  in  the  Avef ria  Dos  calculation  of  Section  3.2, 
it  is  probably  the  case  that  if  turbulence  is  the  cause 
of  structure  it  must  be  localized  in  the  E-layer  coupling 
region. 

Another  consideration  which  must  be  made  in  con¬ 
sidering  neutral  turbulence  as  a  source  of  structure  is 
the  structure  size.  If  structure  were  dependent  upon 
turbulence  in  the  near  altitude  vicinity  of  the  barium 
plasma  then  higher  altitude  clouds  should  either  not 
structure  or  have  larger  structures.  A  dependence  of 
structure  size  upon  altitude  is  not  observed.  However, 
on  the  other  hand,  if  turbulence  is  created  in  the 
E-layer  region  it  definitely  will  control  the  structur¬ 
ing  of  the  barium  cloud  at  higher  altitudes  and  the 
sizes  which  arise  will  be  more  uniform  and  in  closer 
agreement  with  observation,  "ley  questions  of  this 
alternative  explanation  are  is  there  pre-existing 
neutral  turbulence  in  the  E-layer  and  if  there  isn't 
can  the  image  and  image  motion  in  this  region  produce 
neutral  turbulence.  The  investigation  of  these  questions 
is  a  topic  for  future  efforts. 

Clearly  there  are  significant  questions  remaining 
concerning  the  role  of  ion-neutral  coupling  in  the 
morphology  and  structuring  of  electrostatic  ion  clouds. 
These  c[uestions  are  highly  relevant  to  late-time  HANE 
plasma  evolution  because  ion-neutral  coupling  should 
play  a  prominent  role  in  their  evolution.  The  purpose 
of  this  topical  report  is  to  open  up  avenues  of  thought 
into  this  interesting  area  of  HANE  late- time  phenomenology. 
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DEPABTWEHT  OF  THE  AIR  FORCE 
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DEPARTMENT  OF  ENERGY  CONTRACTORS 
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